home *** CD-ROM | disk | FTP | other *** search
/ Sprite 1984 - 1993 / Sprite 1984 - 1993.iso / src / machserver / 1.098 / include / sun3.md / math-68881.h < prev    next >
C/C++ Source or Header  |  1990-09-17  |  9KB  |  390 lines
  1.  
  2. /*******************************************************************
  3. *                                                                  *
  4. *  <math-68881.h>               last modified: March 22, 1989.     *
  5. *                                                                  *
  6. *  Copyright (C) 1989 by Matthew Self.                             *
  7. *  You may freely distribute verbatim copies of this software      *
  8. *  provided that this copyright notice is retained in all copies.  *
  9. *  You may distribute modifications to this software under the     *
  10. *  conditions above if you also clearly note such modifications    *
  11. *  with their author and date.                                     *
  12. *                                                                  *
  13. *  Note:  errno is not set to EDOM when domain errors occur for    *
  14. *  most of these functions.  Rather, it is assumed that the        *
  15. *  68881's OPERR exception will be enabled and handled             *
  16. *  appropriately by the operating system.  Similarly, overflow     *
  17. *  and underflow do not set errno to ERANGE.                       *
  18. *                                                                  *
  19. *  Send bugs to Matthew Self (self@bayes.arc.nasa.gov).            *
  20. *                                                                  *
  21. *******************************************************************/
  22.  
  23. #ifndef __asm
  24. #define __asm asm               /* until new __asm is done */
  25. #endif
  26.  
  27. #ifndef __inline
  28. #define __inline inline         /* until new __inline is done */
  29. #endif
  30.  
  31. #include <errno.h>
  32.  
  33. #ifndef HUGE_VAL
  34. #define HUGE_VAL                                                        \
  35. ({                                                                      \
  36.   double huge_val;                                                      \
  37.                                                                         \
  38.   __asm ("fmove%.d #0x7ff0000000000000,%0"      /* Infinity */          \
  39.          : "=f" (huge_val)                                              \
  40.          : /* no inputs */);                                            \
  41.   huge_val;                                                             \
  42. })
  43. #endif
  44.  
  45. __inline static const double sin (double x)
  46. {
  47.   double value;
  48.  
  49.   __asm ("fsin%.x %1,%0"
  50.          : "=f" (value)
  51.          : "f" (x));
  52.   return value;
  53. }
  54.  
  55. __inline static const double cos (double x)
  56. { 
  57.   double value;
  58.  
  59.   __asm ("fcos%.x %1,%0"
  60.          : "=f" (value)
  61.          : "f" (x));
  62.   return value;
  63. }
  64.  
  65. __inline static const double tan (double x)
  66. {
  67.   double value;
  68.  
  69.   __asm ("ftan%.x %1,%0"
  70.          : "=f" (value)
  71.          : "f" (x));
  72.   return value;
  73. }
  74.  
  75. __inline static const double asin (double x)
  76. {
  77.   double value;
  78.  
  79.   __asm ("fasin%.x %1,%0"
  80.          : "=f" (value)
  81.          : "f" (x));
  82.   return value;
  83. }
  84.  
  85. __inline static const double acos (double x)
  86. {
  87.   double value;
  88.  
  89.   __asm ("facos%.x %1,%0"
  90.          : "=f" (value)
  91.          : "f" (x));
  92.   return value;
  93. }
  94.  
  95. __inline static const double atan (double x)
  96. {
  97.   double value;
  98.  
  99.   __asm ("fatan%.x %1,%0"
  100.          : "=f" (value)
  101.          : "f" (x));
  102.   return value;
  103. }
  104.  
  105. __inline static const double atan2 (double y, double x)
  106. {
  107.   double pi, pi_over_2;
  108.  
  109.   __asm ("fmovecr%.x %#0,%0"            /* extended precision pi */
  110.          : "=f" (pi)
  111.          : /* no inputs */ );
  112.   __asm ("fscale%.b %#-1,%0"            /* no loss of accuracy */
  113.          : "=f" (pi_over_2)
  114.          : "0" (pi));
  115.   if (x > 0)
  116.     {
  117.       if (y > 0)
  118.         {
  119.           if (x > y)
  120.             return atan (y / x);
  121.           else
  122.             return pi_over_2 - atan (x / y);
  123.         }
  124.       else
  125.         {
  126.           if (x > -y)
  127.             return atan (y / x);
  128.           else
  129.             return - pi_over_2 - atan (x / y);
  130.         }
  131.     }
  132.   else
  133.     {
  134.       if (y > 0)
  135.         {
  136.           if (-x > y)
  137.             return pi + atan (y / x);
  138.           else
  139.             return pi_over_2 - atan (x / y);
  140.         }
  141.       else
  142.         {
  143.           if (-x > -y)
  144.             return - pi + atan (y / x);
  145.           else if (y < 0)
  146.             return - pi_over_2 - atan (x / y);
  147.           else
  148.             {
  149.               double value;
  150.  
  151.               errno = EDOM;
  152.               __asm ("fmove%.d %#0rnan,%0"      /* quiet NaN */
  153.                      : "=f" (value)
  154.                      : /* no inputs */);
  155.               return value;
  156.             }
  157.         }
  158.     }
  159. }
  160.  
  161. __inline static const double sinh (double x)
  162. {
  163.   double value;
  164.  
  165.   __asm ("fsinh%.x %1,%0"
  166.          : "=f" (value)
  167.          : "f" (x));
  168.   return value;
  169. }
  170.  
  171. __inline static const double cosh (double x)
  172. {
  173.   double value;
  174.  
  175.   __asm ("fcosh%.x %1,%0"
  176.          : "=f" (value)
  177.          : "f" (x));
  178.   return value;
  179. }
  180.  
  181. __inline static const double tanh (double x)
  182. {
  183.   double value;
  184.  
  185.   __asm ("ftanh%.x %1,%0"
  186.          : "=f" (value)
  187.          : "f" (x));
  188.   return value;
  189. }
  190.  
  191. __inline static const double exp (double x)
  192. {
  193.   double value;
  194.  
  195.   __asm ("fetox%.x %1,%0"
  196.          : "=f" (value)
  197.          : "f" (x));
  198.   return value;
  199. }
  200.  
  201. __inline static const double log (double x)
  202. {
  203.   double value;
  204.  
  205.   __asm ("flogn%.x %1,%0"
  206.          : "=f" (value)
  207.          : "f" (x));
  208.   return value;
  209. }
  210.  
  211. __inline static const double log10 (double x)
  212. {
  213.   double value;
  214.  
  215.   __asm ("flog10%.x %1,%0"
  216.          : "=f" (value)
  217.          : "f" (x));
  218.   return value;
  219. }
  220.  
  221. __inline static const double pow (const double x, const double y)
  222. {
  223.   if (x > 0)
  224.     return exp (y * log (x));
  225.   else if (x == 0)
  226.     {
  227.       if (y > 0)
  228.         return 0.0;
  229.       else
  230.         {
  231.           double value;
  232.  
  233.           errno = EDOM;
  234.           __asm ("fmove%.d %#0rnan,%0"          /* quiet NaN */
  235.                  : "=f" (value)
  236.                  : /* no inputs */);
  237.           return value;
  238.         }
  239.     }
  240.   else
  241.     {
  242.       double temp;
  243.  
  244.       __asm ("fintrz%.x %1,%0"
  245.              : "=f" (temp)                      /* integer-valued float */
  246.              : "f" (y));
  247.       if (y == temp)
  248.         {
  249.           int i = (int) y;
  250.  
  251.           if (i & 1 == 0)                       /* even */
  252.             return exp (y * log (x));
  253.           else
  254.             return - exp (y * log (x));
  255.         }
  256.       else
  257.         {
  258.           double value;
  259.  
  260.           errno = EDOM;
  261.           __asm ("fmove%.d %#0rnan,%0"          /* quiet NaN */
  262.                  : "=f" (value)
  263.                  : /* no inputs */);
  264.           return value;
  265.         }
  266.     }
  267. }
  268.  
  269. __inline static const double sqrt (double x)
  270. {
  271.   double value;
  272.  
  273.   __asm ("fsqrt%.x %1,%0"
  274.          : "=f" (value)
  275.          : "f" (x));
  276.   return value;
  277. }
  278.  
  279. __inline static const double ceil (double x)
  280. {
  281.   int rounding_mode, round_up;
  282.   double value;
  283.  
  284.   __asm volatile ("fmove%.l fpcr,%0"
  285.                   : "=dm" (rounding_mode)
  286.                   : /* no inputs */ );
  287.   round_up = rounding_mode | 0x30;
  288.   __asm volatile ("fmove%.l %0,fpcr"
  289.                   : /* no outputs */
  290.                   : "dmi" (round_up));
  291.   __asm volatile ("fint%.x %1,%0"
  292.                   : "=f" (value)
  293.                   : "f" (x));
  294.   __asm volatile ("fmove%.l %0,fpcr"
  295.                   : /* no outputs */
  296.                   : "dmi" (rounding_mode));
  297.   return value;
  298. }
  299.  
  300. __inline static const double floor (double x)
  301. {
  302.   int rounding_mode, round_down;
  303.   double value;
  304.  
  305.   __asm volatile ("fmove%.l fpcr,%0"
  306.                   : "=dm" (rounding_mode)
  307.                   : /* no inputs */ );
  308.   round_down = (rounding_mode & ~0x10)
  309.                 | 0x20;
  310.   __asm volatile ("fmove%.l %0,fpcr"
  311.                   : /* no outputs */
  312.                   : "dmi" (round_down));
  313.   __asm volatile ("fint%.x %1,%0"
  314.                   : "=f" (value)
  315.                   : "f" (x));
  316.   __asm volatile ("fmove%.l %0,fpcr"
  317.                   : /* no outputs */
  318.                   : "dmi" (rounding_mode));
  319.   return value;
  320. }
  321.  
  322. #ifndef fabs
  323. __inline static const double fabs (double x)
  324. {
  325.   double value;
  326.  
  327.   __asm ("fabs%.x %1,%0"
  328.          : "=f" (value)
  329.          : "f" (x));
  330.   return value;
  331. }
  332. #endif
  333.  
  334. __inline static const double fmod (double x, double y)
  335. {
  336.   double value;
  337.  
  338.   __asm ("fmod%.x %2,%0"
  339.          : "=f" (value)
  340.          : "0" (x),
  341.            "f" (y));
  342.   return value;
  343. }
  344.  
  345. __inline static const double ldexp (double x, int n)
  346. {
  347.   double value;
  348.  
  349.   __asm ("fscale%.l %2,%0"
  350.          : "=f" (value)
  351.          : "0" (x),
  352.            "dmi" (n));
  353.   return value;
  354. }
  355.  
  356. __inline static double frexp (double x, int *exp)
  357. {
  358.   double float_exponent;
  359.   int int_exponent;
  360.   double mantissa;
  361.  
  362.   __asm ("fgetexp%.x %1,%0"
  363.          : "=f" (float_exponent)        /* integer-valued float */
  364.          : "f" (x));
  365.   int_exponent = (int) float_exponent;
  366.   __asm ("fgetmant%.x %1,%0"
  367.          : "=f" (mantissa)              /* 1.0 <= mantissa < 2.0 */
  368.          : "f" (x));
  369.   if (mantissa != 0)
  370.     {
  371.       __asm ("fscale%.b %#-1,%0"
  372.              : "=f" (mantissa)          /* mantissa /= 2.0 */
  373.              : "0" (mantissa));
  374.       int_exponent += 1;
  375.     }
  376.   *exp = int_exponent;
  377.   return mantissa;
  378. }
  379.  
  380. __inline static double modf (double x, double *ip)
  381.                 /* ip means "integer part" */
  382. {
  383.   __asm ("fintrz%.x %1,%0"
  384.          : "=f" (*ip)        /* integer-valued float */
  385.          : "f" (x));
  386.   return x - (*ip);
  387. }
  388.  
  389.  
  390.